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Abstract 

The longitudinal response function R(q, w) of nuclear matter is calculated 
in a semiclassical quark model. The model has a many- body string-flip 
potential that confines quarks into hadrons and avoids color van der Waals 
forces. Molecular dynamics simulations are used to calculate R(q, w) in one 
space dimension for a variety of momentum transfers q and excitation energies 
it;. The response function decreases with density (compared to a free hadron 
response) because of quark exchange between hadrons. 
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I. INTRODUCTION 



The response of nuclei to a variety of electroweak and hadronic probes is rich in excita- 
tions involving many nucleons, such as giant resonances, and in internal excitations of single 
baryons. People have speculated on possible mixing between the "nucleon- like" Gamow- 
Teller resonance and the "quark-like" N — A transition M. We propose to study this mixing 
by calculating a variety of response functions using a quark model of nuclear matter. Indeed, 
the string-flip or quark exchange model of Horowitz et al. and Lenz et al. |J contains 
both types of modes on an equal footing. 

The nuclear excitation spectrum contains many modes with different spin, flavor and 
color quantum numbers. One can speculate on new types of modes involving collective 
excitations of quarks "from" several nucleons. For example, there could be a collective 
oscillation of up quarks against down quarks (independent of weather a given quark is in 
a proton or neutron). This would be qualitatively different from both conventional giant 
resonances (where, for example, neutrons oscillate against protons) and quasifree Delta 
production (which only involves quarks "from" a single hadron). We call such many-quark 
collective modes quark giant resonances and will search for them in our quark model response 
functions. 

We use a string- flip potential model []2|,|3| which (1) confines quarks into hadrons, (2) has 
the correct exchange symmetry among identical quarks (even for quarks "from" different 
hadrons) and (3) has no color Van der Waals forces between color singlet hadrons.0 

The ground state properties of both one dimensional H|J and three-dimensional [|5]|| 
models of nuclear matter have been calculated using variational and path integral (Green's 
function) Monte Carlo techniques. At low density the model reduces to a gas of weakly 
interacting hadrons with properties similar to those in free space, while at high density the 



Such Van der Waals forces are not observed and are present in confining two-body potentials. 
They produce large spurious energies in nuclear matter. 
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model reduces to a free Fermi gas of quarks. 

A first attempt was made [0] to examine the longitudinal response in the model by 
calculating the energy integral of the response function. This integral or Coulomb sum 
can be expressed in terms of the ground state correlation function. The latter was directly 
calculated via Monte Carlo techniques. However, it is difficult to interpret the quark model 
Coulomb sum since it necessarily has contributions from both nucleon excitations (such as 
the quasielastic peak) and quark excitations of a single hadron. To make further progress, 
it is desirable to calculate the distribution of strength in order to separate the nucleon 
excitations at low energy from the "quark" excitations at higher energies. 

Unfortunately, it is not possible to directly calculate the response function in real time 
via path integral Monte Carlo techniques. This is because the complex phases for real time 
evolution introduce prohibitive statistical noise. Some progress has been made calculating 
response functions in imaginary time (see for example ||). However, the imaginary time 
response can not be directly transformed to real time because of statistical noise. Thus, 
there are ambiguities comparing calculations in imaginary time with the physical real time 
response function. 

Alternatively we consider a (semi) classical approximation to our quark model. This 
allows us to calculate the real time evolution directly using molecular dynamics techniques. 
From molecular dynamics trajectories, the response function can be calculated in a straight- 
forward and unambiguous way. 

For our harmonic oscillator string-flip potential (discussed in section II) many of the full 
quantum results are reproduced in a classical approximation. For example, the excitation 
spectrum (multiples of Tiuoq) is unchanged. Furthermore, the ground state correlation func- 
tion g(r) is approximately reproduced. We do this by choosing a fictitious temperature so 
that thermal motion simulates some of the zero point motion of the T = quantum system. 
[Note our calculations are fully classical except for this fake temperature.] 

Our classical approximation may allow insight into any collective modes we observe. 
New modes with a classical analog may be simple to understand. We expect our classical 
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approximation to be better at describing long wavelength modes such as zero sound rather 
than short wavelength modes such as the quasielastic peak. Therefore, we will focus on the 
long wavelength response. 

In this first study we consider a very simple model in one spatial dimension without 
explicit spin, color or flavor degrees of freedom. This same model was considered in Ref. 0. 
It is straightforward to extend our calculations to three dimensions (see for example f|). 
Furthermore, one can include the internal degrees of freedom and add additional two-body 
interactions in order to reproduce the two nucleon phase shifts |§. We hope to extend our 
simple model in later work. 

Our model is described in section II along with the molecular dynamics formalism. Sec- 
tion III describes the free quark and free hadron response functions. Next, results for nuclear 
response functions at various densities, momentum transfers and excitation energies are pre- 
sented in section IV. Finally, we summarize in Section V. 

II. FORMALISM 

In this section we describe the string-flip potential model and then the molecular dy- 
namics formalism for calculating the response function. Our model has two-quark hadrons. 
Although simple, these have interesting composite features. It is straightforward to extend 
our calculations to three quark hadrons ||. 

Our classical approximation does not enforce antisymmetry in the spatial quark wave 
function. However, the relatively large spin-flavor and color degeneracy of nuclear matter 
suggests that most of the time the quarks are antisymmetrized in the internal coordinates. 
Therefore, the wave function is symmetric in the spatial coordinates. We hope to improve 
the treatment of internal degrees of freedom in future work. For now, we consider the 
simplest model without explicit spin, flavor, or color degrees of freedom. 
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A. String Flip Model 



The Hamiltonian of our model is given by (for a one dimensional system) , 

H = T + V, (1) 

where, 

-k 2 N d 2 

N/2 

V = min {P} v{x P (2n-i) - x P{2n) ). (3) 

n=l 

Here the minimum is taken over all permutations P of the quark labels which pair the N 
quarks into N/2 pairs, with the confining potential v or link acting only between the two 
quarks in each of the N/2 pairs. This insures that all quarks will be confined into hadrons 
and that there will be no color Van der Waals forces between hadrons. Furthermore, the 
minimum prescription insures the potential is symmetric in all quarks coordinates. For the 
four quark system (for example), the potential is, 

V(x 1 ,x 2 , x 3 , x 4 ) = minjVi, V 2 , V 3 }, (4) 



V\ = v(x 12 ) + v(x 34: ), 



V 2 = v(xi 3 ) + v(x 24 ), 



(5) 



(6) 



V 3 = v(x 14 ) + v(x 23 ), (7) 

Thus, although the potential is an N-body interaction (i.e., it depends upon the configuration 
of all N quarks), it is basically quite simple. The confining forces operate only between pairs 
of quarks. At any instant, the N/2 links arrange themselves in a way that minimizes the 
potential energy; this idea is familiar in adiabatic approximations to field theories. The only 
hadron interactions come from the exchange of quarks between clusters. 
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The minimum pairing for a large N system (in three space dimensions) can be determined 
with special techniques (see for example || and |I0|). However, in one dimension (assuming 
periodic boundary conditions) this problem is simple. Simply order the coordinates from 
left to right and then the minimum is one of only two possibilities, 

V = mm(V R ,V L ), (8) 

where 

V R = \Y, (^~^+i) 2 , (9) 

i=odd 

V l = \Y. (x i -x i+1 ) 2 + -{x N -xx -Lf. (10) 

i=even ~ 

Here L represents the length of the box. Any other pairing will involve overlapping links 
and have higher energy. Note, for simplicity we have assumed harmonic links, 

v(r) = (^)r 2 . (11) 

We discuss below our choice for the oscillator constant k. 

This model Hamiltonian has been used in previous works for quantum mechanical cal- 
culations The expectation value of the energy E is calculated for a trial wave function 
(note as discussed before, we are assuming Bose symmetry for the spatial coordinates), 

V=exp{-\V}. (12) 

Here V is the full many-body potential and A is a variational parameter. The value of 
A is chosen such that E = (ty\H\ty) is a minimum for a given density p of our system. 
At low densities, A = l/v2 reproduces the oscillator wave-function for isolated hadrons 
while at high densities A — > reproduces a free particle system. It is found that as the 
density increases, the value of A that minimizes the energy decreases. The value of A for 
each density will be used to determine the temperature of our many-body system for the 
molecular dynamic simulation (see below). 

6 



Our results can be readily scaled for different values of m and k. For convenience we 
chose the length and energy scales such that the quark mass to, the spring constant k and 
the Planck's constant H are all equal to one. The energy scale is determined by the energy 
of a free hadron and is given by 

E = (n + ^)htu, (13) 

with 

hu = hJ—. (14) 

V to 

The length scale is determined by the root mean square seperation of the quarks in the 
ground state hadron and is given by, 

/oo 
p{x)x 2 dx, (15) 
-oo 



where 



P( " ) = ( 2^ )V4eXp{ - ( ^ )1/V} - (16) 
After evaluating the integral we get, 

<* 2 > 1/2 ^> 1/4 = omi O 1 "- < i7 > 

These are our scales for energy and length in our calculation. They might be set to 
describe the radius of a nucleon, or a fraction of the charge radius described by a quark core 
(neglecting the pion cloud) and a typical baryon excitation scale. 

We determine A by minimizing the expectation value of the energy E. Note that 

E = (2\ 2 + 1)(V) X: (18) 

where 

JVe~ 2XV dT 

oo. = TT^r- (19) 



The equation for E was obtained by using Eqs. (1,2,8-10,12). In principle, one could 
proceed and calculate 4^ , setting the resulting expression equal to zero and determining A 
numerically such that the derivative is zero. But there is a simpler way of doing this using 
the following scaling property 0: 

(V) P ,x = (^) 2 (V) p ,,y, (20) 

This property is used to get expressions for A' , p' and E(p') in terms of two input parameters 
which we denote by A and p. 

Summarizing, the prescription is as follows : give A and p as two input values. Then 
calculate A', p' and E(p') using the equations given in ref ||. The value of A that minimizes 
the energy E is our A'. 



B. Molecular Dynamics Simulation 

A molecular dynamics simulation determines the motions of the particles in a system 
PH . It is based on the microcanonical ensemble and assumes ergodicity. Our many-body 
potential, eq. (3) or (8), is used in conjunction with Newton's equations of motion to 
determine the trajectories of the particles. 

The steps of our molecular dynamics simulation are the following : 
a) The initial coordinates of the particles are determined by using the partition function 
in a Metropolis simulation. The partition function is given by, 

Z = exp{-V/k B T}, (21) 

where k B is the Boltzmann constant and T is the temperature of the system. The temper- 
ature of the system is determined by the prescription, 

k B T = JL (22) 

where A was determined in a quantum variational calculation (see fig. la). This prescription 
for choosing a fake temperature insures that the classical calculation will reproduce the 
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correlation function g(r) and other ground state observables of the variational quantum 
mechanical calculation. That is, instead of sampling the square of the wave function eq. 
(12) we sample an equivalent partition function eq. (21) in the classical calculation. This 
partition function is used in a succession of Metropolis sweeps to evolve the system until it 
reaches equilibrium . We used 1000 Metropolis sweeps to determine the initial coordinates. 
This procedure simulates quantum zero point motion with thermal motion. 

b) The initial velocities of the particles are randomly selected from the Boltzmann dis- 
tribution given by, 



l2nk B T -\mv 2 . 

c) The trajectories of the particles are determined by numerically integrating Newton's 
equations of motion for which we used an algorithm by Beeman [pl|. Here a time step of 
0.01 insures that energy is conserved to better than one percent (see fig. lc). 

d) From the trajectories we can calculate any time-dependent quantity of interest. The 
time-dependent quantities considered here are the following : 

(i) One-body density, 

1 N 

Q(q,t) = ^=J2 eiqXj(t) ( 24 ) 



This represents the Fourier transform of the charge distribution of the system where all the 
quarks are assumed to be point-like, with unit charge. In a later work we will consider 
quarks of different flavors and charges. 

(ii) Two-body density or radial distribution function, 

(<7(r-))time, (25) 

where 

9(r)= ^^^ E^-I^-^D- ( 26 ) 

Here p is the density of the system, N is the number of particles and T{ and rj are the 
ith and jth particle positions , respectively. This correlation function is proportional to the 
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probability of finding a particle at a distance r from another given particle and is normalized 
to one for large r. 

(iii) Autocorrelation function, 

F(q,t) = (- [ T Q(q,s) Q*(q,s + t) ds). (27) 

T JO 

The brackets mean that we take an ensemble average, that is, an average over many different 
initial conditions. For t — F(q, 0) = S(q) where S(q) is the static structure factor. S(q) 
can also be expressed as, 

N r 

S(q) = 1 + — / dr e iqr [g(r) - 1]. (28) 
From eq. (27), one can show that, 

F(q,-t)=F*(q,t), (29) 

which will be useful in simplifying the calculation for the response function. 

(iv) Response function or dynamical structure factor, 

1 roc 

R(q,w) = — dte lwt F(q,t). (30) 

Z7T J~oo 

By using the property F(q, —t) = F*(q,t) we can rewrite R(q,w) as, 



R(q,w) = -Re dte lwt F(q,t). (31) 

7T JO 

Thus R(q,w) is real. Of course the integral of R(q,w) gives S(q), 

/oo 
dw R(q,w). (32) 
-oo 

The response function is the goal of our calculation since it will provide us with information 
about the excited states and collective modes of the many quark system. 

C. Autocorrelation formalism 



The longitudinal response function R(q, w) describes the response of a system to a weak 
external probe. In inelastic electron scattering R(q, w) is proportional to the differential 
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cross section da ||14|| . An electron probe will interact with the charge density of the target 
nucleus (we ignore the transverse response for now). Therefore, the interaction Hamiltonian, 
Hi n t, will be proportional to the charge density, p(x). The charge density is given by, 



N 



P( X ) = e ^ ( X ~ X i) ( 33 ) 

i 

where for our calculation e» = 1. One can show (by using Fermi's golden rule) that da oc 
R(q, w) where, 

R(q,w) = Y,\(f\pM 2 5(E f -Ei-w), (34) 
/ 

is the quantum mechanical definition for the response function. Here p q = J2j eje tqXj = Q is 
the Fourier transform of the charge density p(x). 

The quantum mechanical R(q, w) can be rewritten as, 

/OO (If 
_^ wt (?\p\{t)p q {m- (35) 

Now, in our case we are dealing with a classical system, so R(q, w) becomes, 

/°° d,t 1 f T 

-e^(- dt'p* q (t + t')p q (t')). (36) 
-oo Z7T T JO 

Thus we are calculating in eqs. (27,31) the classical analog of the quantum mechanical 
response function. 

III. FREE QUARK AND FREE HADRON RESULTS 

We will compare our nuclear matter results to the response of a gas of free quarks and 
to the response of a collection of free hadrons. Therefore in this section we derive analytic 
expressions for these two responses. 

A. Free Quark Response Function 



The trajectories of a free quark system with no interactions are straight lines and con- 
sequently we can find a closed expression for the response function. Suppose we have N 
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quarks, where each of them has a given initial velocity v Q i and a given initial position x {. 
Then, the trajectories are, 

Xi(t) = x oi + v oi t. (37) 
This implies that Xi(s +t) — Xi(s) + v Q it which then gives us, using eq. (24), 

Q(q,s)Q*(q,s + t) = ±J2Y, eiqXi{s)e ~ iqXj{a+t) - (38) 

» 3 

Now, when we integrate over the time variable s in order to calculate the autocorrelation 
function, only the case i — j gives a contribution. Therefore, 

n<i,t) = jjTL*-* *- (39) 

j 

By using our expression for R(q,w) we get that, 

R(q,w) = j-J25(w-qv 3 ). (40) 

j 

We take the average of this expression for a Boltzmann velocity distribution to give, 

R(q,w) = r dvS(w - qv) J_ e- 2 /^ T ). (41) 

J-oo \J11TKb1 

Performing the integration we finally get [I2~|, 



R(q, w) = - \= e - w2 'W kBT \ (42) 

This is our expression for the response function of a free quark system. To get the autocor- 
relation function F(q,t) simply Fourier transform eq. (42), 

F(g,t)=exp{-^^}. (43) 



B. Free Hadron Response Function 

The free hadron case implies that two quarks are interacting via a harmonic oscillator 
potential and there is no exchange of partners. For this case we know the exact trajectories 
for each of the two quarks and can calculate the response as shown in the appendix. 
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IV. RESULTS 



We now present and discuss the results obtained in our calculation. In Fig lb we plot 
sample trajectories for N = 8 quarks as a function of time for a density of p = 0.414. Note 
that periodic boundary conditions are used with a box length of L = 19.32 (which is N/p). 
This figure illustrates many of the features of the model. 

At t=0 the 8 quarks are grouped into four separate hadrons by the string-flip potential. 
For an isolated hadron, the potential is simply a harmonic oscillator, so that the trajectories 
consist of two curved lines oscillating together. The average slope of the two lines gives the 
center of mass velocity of the hadron and the amplitude of the oscillations gives the internal 
excitation energy. Near x=0 and t=32 we see a hadron move through the left wall of the 
box. Since we are using periodic boundary conditions the two quarks reappear near x=19 
at the far wall of the box. 

A variety of different kinds of hadron-hadron collisions are possible with our potential 
via string rearrangement. For example, near x=5 and t=23 there is an inelastic collision 
between a rapidly downward moving hadron and a slower upward moving hadron. After 
the collision, the hadron which emerges at larger x has a larger internal excitation energy 
which is indicated by the larger oscillation amplitude. Clearly this energy came from the 
center of mass motion of the two original hadrons. Note that hadrons only interact via 
string rearrangement. If two quarks are very close together, it is likely they will remain 
paired so that the probability to interact is small. For example, near x=5 and t=87 there is 
a hadron moving upwards with the two quarks very close together. As result the quarks are 
only slightly perturbed when this hadron passes through a larger downward moving hadron. 
This behavior follows from our assumption about the saturating nature of the forces (once 
two quarks are paired off the remaining interactions are zero). It may be the analog, in our 
very simple model world, of QCD inspired ideas regarding color transparency. 

In Fig.lc we also show the total energy of the N=8 quark system vs. time. Our inte- 
gration method with a step size of dt=0.01 conserves energy to better than one percent. 
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Furthermore, there is little evidence of a long term gain or loss in energy. The hardest part 
of determining the trajectories is finding when the strings flip. We simply check for flipping 
at each time step. No attempt has been made to try and evolve the harmonic oscillators 
analytically between string-flips. 

In fig. (2) we show the radial distribution function g{r) calculated in two completely 
different ways. The large value of g(r) for small r represents the other quark which is 
bound into a given hadron by the string-flip interaction. The dots refer to a Metropolis 
calculation which could either represent a quantum Monte Carlo sampling of the variational 
wave function of eq. (12) or a sampling of the classical partition function of eq. (21). The 
stars refer to a Molecular Dynamics simulation. We get good agreement between the two 
methods. This provides a check of our Molecular Dynamics numerics and illustrates our 
simulation of quantum zero point motion with a fictitious temperature. In Fig. (3) we 
show the autocorrelation function F(q, t) as a function of time. The curve with the stars 
corresponds to a free hadron calculation, the curve with the dots corresponds to the string- 
flip potential and the solid curve corresponds to the free quark case. We recall that from 
the definition of F(q,t), this function measures how the value of Q at s + t is correlated 
with its value at t. If particles are moving very slowly (low velocities), then F(q,t) will 
fall off, as t — > oo, very slowly since the value of Q at s + t is very correlated to its value 
at t. The opposite will happen if the particles are moving very rapidly. We see that the 
free quark F(q, t) falls off faster than the other two autocorrelation functions. The single 
hadron F(q, t) is generally above the corresponding curves for the free quark and the string 
flip potential. We note that the string-flip F(q, t) oscillates slightly negative before going to 
zero at large t. The response function R(q, w) is shown in Fig. (4) for a low density. At this 
low density, the string-flip response function is very close to the response for free hadrons. 
As the density is increased in Figs. (5-8) we see that the string-flip response becomes smaller 
than the free hadron result (for low excitation energy w) and eventually approaches the free 
quark response at high densities. 

This decrease in the string-flip response with density illustrates the decrease in the co- 
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herence of scattering from the two quarks in a hadron (without break-up). At low density 
and momentum transfer one scatters coherently from both quarks. However, as the density 
increases, it becomes possible for quarks to hop from one hadron to another because of the 
string-flip potential. This suggests that quark exchange effects will decrease the response of 
nuclear matter at low excitation energies. We will discuss this more in the next section. 

In Fig. (6) we show the response function R(q, w) for different momentum transfers q. 
Results at high q look qualitatively similar to the low q responses. This may be an artifact 
of our semiclassical approximation. We expect the classical results to be better at low q. 

Finally, in Figs (9-10) we test for finite size effects by comparing response functions 
calculated with N=8 and N=16 quarks. At a low density of p = 0.211, finite size effects are 
very small. They are somewhat larger at p = 0.414, but they are still small. Therefore we 
do not expect finite size effects to change our N=8 results greatly. 

V. SUMMARY AND CONCLUSIONS 

In this paper we have presented a semiclassical simulation of a simple quark model of 
nuclear matter. We have used thermal motions from a fictitious heat bath to simulate quan- 
tum zero point motion. Our string-flip potential model: (1) confines quarks into hadrons, 
(2) allows the hadrons to separate without color Van der Waals forces and (3) is symmetric 
in all of the quark coordinates. The low density limit of the model is a gas of free hadrons, 
while the high density limit is a free quark gas. This string-flip model has been used previ- 
ously for a variety of quantum calculations of nuclear matter ground state properties. We 
have calculated the longitudinal response function from a molecular dynamics simulation. 
We find that R(q, w) progresses smoothly from that for free hadrons to free quarks as the 
density increases. Compared to the response of free hadrons, R(q, w) is found to decrease 
with increasing density (at low w) . This represents a decrease in coherence as quarks "hop" 
from one hadron to the next, i.e., at low momentum transfers one scatters coherently from 
all of the quarks in a free hadron. However, the possibility of quark exchange between 
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nucleons decreases this coherence in nuclear matter. Our results suggest that this quark 
exchange effect may somewhat decrease the strength of traditional hadronic calculations of 
the nuclear response. 

Indeed, people have speculated that quark effects may decrease the Gamow Teller 
strength, while the experimental integral of the longitudinal response (the so called Coulomb 
sum) could also be below theoretical expectations. It would be interesting to explore this 
suggestion in more sophisticated quark models. 

Our simulation provides no evidence for possible new collective modes. However, the 
simulation has been for the simplest model without explicit spin, flavor, or color internal 
degrees of freedom. In future work, we will include more internal degrees of freedom (see 
for example ||15|| ) and this will allow us to search a much richer excitation spectrum for new 
collective modes with various spin, flavor, and color quantum numbers. 

For simplicity, our calculations have been in one space dimension. However, three di- 
mensional simulations should be straightforward (although requiring more computer time). 
Efficient algorithms do exist for evaluating the string-flip potential in three dimensions PJIU|. 
Finally, we plan to carry out quantum mechanical calculations of the response function in 
imaginary time via Green's function Monte Carlo to compare with these semiclassical results. 
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In this appendix we calculate the response function for two isolated quarks interacting 
via a harmonic oscilator potental Let the coordinates be X\ and Xi- The potential energy 
between the two quarks is given by, 



APPENDIX A: FREE HADRON RESPONSE FUNCTION 




(Al) 
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Introducing a relative coordinate x and a center mass coordinate x cm we can write x\ and 

x 2 as, 

1 

X\ — x cm ~ ^ (A2) 
X 2 = Xcm + ^X. (A3) 

The trajectories x cm (t) and x(t) are given by, 

x C m(t) = A cm + -B cm t (A4) 
x(t) = A cos(v / 2t + 5). (A5) 

where A cm , B cm , A and 5 are constants. Substituting these results in the expression for 
F(q,t) and after some simplification and a change of variable we get, 

F(q, t) = -(e- 1 ^ 1 cos[A^cos(y)]cos[A^-cos(y + V2t)]dy). (A6) 

7T JO 2 2 

We can simplify this even more by defining the function g(t) as, 

g(t) = f W cos\A^cos(y)]cos[AUos{y + V2t)]dy, (A7) 
JO z z 

were g(t) is periodic with period ir/ Consequently we can expand g(t) in a Fourier series 
as follows: let — oo < t < oo with L = , then g(t + 2L) = g(t) is 2L periodic and 



therefore, 



g(t) = A + jr(A n cos(^) + B n sin(^)). (A8) 
n=i L L 

Note B n = because g(t) is an even function in t, that is, t) = g(t),so that 

g(t) = A + J2A n cos{—). (A9) 

n=l L 

The coefficients A and A n are given by, 

1 r L 
LJo 

2 r L _ ,nvrt 



^o = T / <7(*)<ft, (A10) 
L Jo 

2 f L , x ,n7rt. , . , „ . 

n= LJo 9(t)cos(—)dt, (All) 
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Inserting these results in the expression for F(q,t) given above we get that, 

O O CO 

F(q,t) = -(e- tqVcmt A ) + _( e -^-* V A n cos(2V2nt)) . (A12) 

n * n=l 

Recalling that R(q, w) is the Fourier transform of F(q, t) we get after substituting eq. (A12) 
into eq. (30), 

2 

R(q,w) = -(5{w - qv cm )A ) + 

7T 

j oo 

- ]T{(,4 n <5(u> + 2^2n - gw cm )) + (A„<5(u> - 2y/2n - qv cm ))}. (A13) 

71 n=l 

Let's write explicitly what the term in brackets means, so that we can see how this 
expression is calculated, 

The equations for (A n 5(w ± 2^^ — qv cm )) are of the same form. Finally, by simplifying 
and making a change of variables in the integrals that represent the expressions in brackets 
we get the following expression for the response function: 

2 e -w 2 /(q 2 k B T) fn / 2 



R(q,w) = — ti = w~ \ dz d y 

n 2 q J-KknT Jo Jo 



ir 2 q V^keT 

k T k T 
x ( ex P[ |-9 2 (cos(y) - cos(y + z)) 2 } + exp[ ^q 2 (cos(y) + cos(y + z)) 2 ]} 

i - 2 v { [ e -(^)V(^T) + e -(^)V( fcB T) ] r' 2 dz r dy 

n 2 q\fnk B T Jo Jo 

x ( ex P[ |-9 2 (cos(y) - cos(y + z)) 2 } + exp[ 1— g 2 (cos(y) + cos(y + z)) 2 ]) 

x cos(2nz)}. (A15) 

This is our final expression for the response function of a single two-quark hadron. It looks 
complicated, but it is straightforward to calculate numerically since we have a double integral 
of a smooth function. 
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FIGURES 



la) Value of the variational parameter A which minimizes the energy for a given density 
p. Results are for a quantum variational Monte Carlo calculation involving 8 quarks using 
the wave function in Eq. (12). 

lb) Trajectories for 8 quarks. 

lc) Energy conservation in the algorithm used for the integration of Newton's equation 
of motion. This is the total energy for 8 quarks as function of time. 

2. Radial distribution function for p = 0.414 calculated in two different ways. The dots 
and the stars represent a Metropolis and a Molecular Dynamics calculation, respectively. 

3. The autocorrelation function versus time. The curve with the stars corresponds to 
the free hadron case, the curve with the dots and the error bars corresponds to the string 
flip model and the solid curve correspond to the free quark case. The values of momentum 
transfer and density are indicated at the top of the graph. 

4. The response function R(q, w) as function of frequency w for a momentum transfer 
of q — 0.034 and a density of p = 0.043 . The curve with the stars corresponds to the 
free hadron calculation, the curve with dots and error bars corresponds to the string flip 
potential and the solid curve corresponds to the free quark case. 

5. The same as in Fig. 4 but with q = 0.33 and p = 0.212 . 

6. The response function R(q,w) with a fixed value of density equal to p = 0.414 and 
three different values for the momentum transfer q . a) The same description as in Fig. 4 
but with q = 0.32 and p = 0.414 . b) The same description as in Fig.4 but with q = 0.65 
and p = 0.414 . c) The same description as in Fig.4 but with q = 1.30 and p = 0.414 . 



21 



7. The same as in Fig.4 but with q = 0.70 and p = 0.893 



8. The same as in Fig.4 but with q = 1.25 and p = 1.59 . 



9. Finite size effects in the response function R(q,w) for q = 0.17 and p = 0.21 . The 
curve with the dots corresponds to a calculation with 8 quarks and the curve with the stars 
corresponds to a calculation with 16 quarks. 



10. The same as in Fig.9 but with q = 0.32 and p = 0.414. 
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